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in 

We numerically study the binary disordered Poland-Scheraga model of DNA denaturation, in the 
' regime where the pure model displays a first order transition (loop exponent c = 2.15 > 2). We 

CI use a Fixman-Freire scheme for the entropy of loops and consider chain length up to N = 4 • 10 , 

with averages over 10 4 samples. We present in parallel the results of various observables for two 
boundary conditions, namely bound-bound (bb) and bound- unbound (bu), because they present 
very different finite-size behaviors, both in the pure case and in the disordered case. Our main 
conclusion is that the transition remains first order in the disordered case: in the (bu) case, the 
ly-^ ' disorder averaged energy and contact densities present crossings for different values of N without 

^s^j , rescaling. In addition, we obtain that these disorder averaged observables do not satisfy finite size 

scaling, as a consequence of strong sample to sample fluctuations of the pseudo-critical temperature. 
For a given sample, we propose a procedure to identify its pseudo-critical temperature, and show 
that this sample then obeys first order transition finite size scaling behavior. Finally, we obtain that 
the disorder averaged critical loop distribution is still governed by P(l) ~ l/l c in the regime I <C N, 
iy>j ' as in the pure case. 



£ 

I 



O 

o 



(N 



O 

c3 



X 



I. INTRODUCTION 



Recent Monte Carlo simulations of three dimensional interacting self avoiding walks (SAW's) Ql H, Q| have revived 
interest [HIS 00 m the Poland-Scheraga (PS) model of DNA denaturation [jj. In this model, the central parameter 
is the exponent c governing the weight \/l c of a loop of length I. If the binding energy is constant along the chain, 
the model is exactly soluble, and the transition is first order for c > 2 and second order for 1 < c < 2. Gaussian loops 



in d dimensions are characterized by c = cq = i. The role of self avoidance within a loop was taken into account 



by Fisher and yields c = dv : where v is the SAW radius of gyration exponent. More recently, Kafri et al. 0,0 
pointed out that the inclusion of the self avoidance of the loop with the rest of the chains further increased c to a 
value c > 2, both in d = 2 and d — 3 (see also ^l]). The transition in the pure model should then be discontinuous, 
as previously found by in Monte Carlo simulations of SAW 's. The value c ~ 2.11 was in turn measured in three 
<^ \ dimensional Monte Carlo simulations by 

. The precise correspondence between the PS model and the numerical simulations of SAW's in d dimensional space 
is a matter of debate: the comparison of the finite-size properties of both models |8| involves not only the loop weight 
l/l c discussed above, but also the end segment properties, since simulations use SAW's bound at the origin, but 
with free end points. The PS model can easily be formulated for free boundary conditions, but this requires some 
prescription for the weight of the unbound free ends. The correct prescription to mimic the SAWs model is still 
controversial |3j,|6|,|8j as explained in more details below. 

In this paper, we are interested in the effect of disorder on the PS model with loop exponent c = 2.15. This disorder 
approach is a first step towards the heterogenity of biological sequences, stemming from the existence of AT and GC 
Watson-Crick pairs. We will therefore consider a binary disorder distribution of the binding energy. In a physics 
oriented context, there have been few studies of the effect of disorder on PS and PS-related models with c > 2. For 
SAW's, Monte-Carlo simulations have been interpreted as follows (i) in 0, the loop exponent c > 2 was found to 
remain the same as in the pure case for two particular biological sequences (ii) the more complete study of the binary 
disordered case in ref ^2] suggests a smoothening out of the transition, without fully excluding the possibility of a 
first order transition. 

On the theoretical side, there exists many discussions on the effect of disorder on first-order transitions fl^lF4llT5llTr3 | 
but these studies typically consider spin systems displaying coexisting domains in the pure case, so that the conclusions 
of these studies cannot be directly applied to a polymer transition, such as the present PS transition. In this paper, we 
study numerically the disordered PS model for a binary distribution of the binding energies. We use the same method 
as in our recent study |17| of the disordered two-dimensional wetting model (corresponding to the PS model with 
c = 3/2). We present in parallel results for two different boundary conditions bound-bound (bb) and bound-unbound 
(bu), since these two boundary conditions lead to very different finite size behaviors, as already stressed in 8] for the 
pure case. 
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II. POLAND-SCHERAGA MODEL OF DNA DENATURATION 



A. Model and observables 



We consider here pairing interactions between monomers on the two chains, in contrast to the stacking interactions 
present in the usual treatment of DNA melting 0, 0, [2(j • The Poland-Scheraga description of DNA denaturation [t| 
can be defined in the following way: we consider a forward partition function Zf(a) for two chains of (a) monomers, 
with monomer pairs (1) and (a) bound together. 




FIG. 1: Recursion relation for Zf(a) (eq.Q) in the PS model with pairing interactions e Q 
This partition function obeys the recursion relation (see Fig^) 



Z f (a)=e-^ J2 Z f (a'W(a';a) (1) 

a' = l,3,.. 

where A/"(a'; a) is the the partition function of a loop going from a' to a. In polymer physics, the asymptotic expansion 
of Af(a'; a) for (a — a' 3> 1) is given by |21| 

A/V; ") - fJ> a a ' f(a - a') (2) 

where <tq is a constant, fi is the fugacity of a pair of monomers and f(x) = i is the weight of a loop of length x ^> 1. 
In the PS model, equation is assumed to hold at all loop lengths. The exponent c is an input of the model, and 
we study here the case c = 2.15, as explained in the Introduction. In a similar way, we define a backward partition 
function Z b (a), defined as the partition function of the two chains starting at monomers (TV), with monomers (a) 
bound together (monomers (TV) being bound or unbound). It obeys the recursion relation 

Z b (a) = e-P e " I Z b {a')M{a- a') + M(a; N) ) (3) 

where 

M(a;N) ~ oi g(a;N) (4) 

is the weight for the possible unpaired segments of length (N — a), and where <r\ is a constant. The first term on 
the r.h.s of equation J2J is similar to the one of equation The second one allows for the boundary condition at 
monomer (iV). For (bb) boundary conditions, a% = 0, since monomers (N) are bound together. For (bu) boundary 
conditions, o\ ^ 0, and the function g(a;N) is an additional input of the model. According to Kafri et al. [g, one 
has g(a; N) ~ ( . jv _ 1 a ^ ei for (TV — a <§; TV), with c\ ~ 0.092 for d — 3. The precise form of g(a; TV) beyond the regime 

(TV — a <§; TV) is still controversial 0,13- Our choice is this paper is to take g(a;N) ~ ny^tpi anc ^ c i = f° r au 
values of a, as is done in biological use of the PS model. Other values of ci and/or g(a; TV) will be briefly considered 
in the pure case. 

In these notations, the thermodynamic partition function Z is given by Z = Z b (l), and the probability for monomers 
(a) to be bound is 

= Z f (a)Z b (a)e^ 
Z b (t) 
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where the factor e^ e ° in the numerator avoids double counting of the contact energy at a. We will also use the 
impairing probability q{a) — 1 —p{a). A quantity of primary importance in the DNA context is the fraction of paired 
monomers, or contact density 



1 N 



a=l 



In the pure case, 9n{T) is proportional to the energy. Since this is not true in the disordered case, we also consider 
the (contact) energy density 



cn{T) = ^-^Te a p(a) (7) 



a=l 



We will also be interested in Pi oop (a, 7), defined as the probability of having a loop between bounded monomers at 
a and 7 



p ( n <v\ ^/(Q)A^(a; 7)^(7) . s , 

-floop [a, 7) = ~z^X) ' 



As in (Tj|, we define the probability measure Mjsr(l) for the loops existing in a sample of length N as follows: for 
each loop length I, we sum over all possible positions (a, 7 = a + I) of the loop probability Pi 00 p(o!, a + I) defined in 
(eq.@) 



N-l 



(a, a + l) (9) 

a=l 

The normalization of this measure over I corresponds to the averaged number of loops in a sample of size N, or 
equivalently to the averaged number N9n(T) of contacts © between the strands. 



B. Numerical implementation 



The above equations show that the numerical calculation of the thermodynamic partition function Z = Zb(l) 
requires a CPU time of order 0(N 2 ). The Fixman-Freire method reduces this CPU time to O(N) by approximating 
the probability factor f(x) of equation by 

1 1 

f(x) = — ~f FF (x)=Y. a >< e ~ bkX for 1<^<^ (10) 
fe=l 

This method is very well known and widely used in biology, since the standard program MELTSIM, which yields 
melting curves of DNA sequences, uses the Fixman-Freire scheme with / = 14 terms. This procedure has been tested 
on DNA chains of length up to N = 10 6 base pairs [TIE^. 

Since the replacement (|1(J|) can be rather surprising for physicists, we "justify" it in a more detailed way. One point 
of view is to start with the integral representation of the T function 

1 1 r M 

- = liniM^— - / dt f - 1 e~ u (11) 

i° r ( c ) Jo 

and to discretize it for L — Ml large as 

^ m 

— = lim ro ^oo a k e ~ bkL ( 12 ) 

This shows that a power-law can be represented by a discrete sum of exponentials if their number m, positions bk and 
weights Ofc are conveniently chosen. In the Fixman-Freire procedure, the 21 coefficients (o^, b^) are obtained from the 
fit on 21 points (Zj., .., Zaz) such that (lnZ 1 ,lnZ 2 , •■,ln/2/) divide into equal intervals the domain [lnZ m j n , hil max ] where 
Irain and Imax are respectively the minimal and maximal loop lengths that are needed in the numerical program. 
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TABLE I: The coefficients list of the Fixman-Freire scheme (eq. 1101 1 used in this paper. 



The number 21 of coefficients is then chosen to obtain the desired numerical accuracy. It turns out that the choice 
7=15 gives an accuracy better than 0.3%. We have adopted this value throughout this paper, with the values of the 
coefficients (a,k,bk) given in Table|U The reader is invited to draw the log-log plot of the curve /ff(0 = a fe e_bfc/ 
from l m in = 1 up to l max = 8 ■ 10 5 (maximum value used in this paper): it turns out that is indistinguishable from a 
straight line with slope —2.15. 

Putting everything together, the model we have numerically studied is defined by recursion equations (lll.'il) for 
the partition functions where (i) the loop partition function AT (a; a') has been replaced by its asymptotic expression 
(J2J, with the Fixman-Freire approximation (eq. 1 1 (J|) for f(x) (ii) the end segment partition function M.{a;N) has 
been replaced by its asymptotic expression (@J, with g(a; N) = g(N — a) = ny-Hpi > w i tn c i = 0- We have taken 
do = 0.29607 and a\ — 0.5 for (bu) boundary conditions. 

Before we turn to the disordered case, we first discuss the finite-size properties of the pure model that will be useful 
to interpret the results of the disordered case. 



III. FINITE SIZE PROPERTIES OF THE PURE CASE 



A. Contact density 




FIG. 2: The fraction of paired monomers 8n(T) of the pure case for sizes N = 
The boundary conditions are (a) bound-bound (b) bound-unbound. 



10 5 (O), 2 ■ 10 5 (□), 4 • 10 5 (o), 8 • 10 5 (A). 
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FIG. 3: The scaled fraction of paired monomers 6((T - T C )N) of the pure case, for T c = 64.521.. for sizes N = 10 5 (O), 2 ■ 10 5 
(□), 4 ■ 10 5 (o), 8 ■ 10 5 (A). The boundary conditions are (a) bound-bound (b) bound-unbound. 



We choose e a = So = —355 K (see below). The contact density 0n(T) defined in cq. J()J) is plotted on Fig. [21 for 
the two boundary conditions (bound-bound) and (bound-unbound) : for the (bb) case, there is no crossing at T c as 
N varies, whereas for the (bu) case, there is a crossing at T c ~ 64.521.. (Celsius scale), as was also found for the pure 
SAWs Monte-Carlo simulations (see Fig. 7 of Ref flj])- From a numerical point of view, the boundary conditions (bu) 
seem therefore much more interesting since it yields a direct measure of the critical temperature T c . With the value 
of T c obtained from the crossing of Fig. 1 (b), we may now plot on Fig. 2 the contact density 6n{T) in terms of the 
rescaled variable (T — T c )iV* with the crossover exponent cj) = 1 : for the two boundary conditions (bb) and (bu), the 
finite-size scaling form 

9 N (T) = e((T-T c )N) (13) 

holds in the delocalized phase T > T c , as was also found in the pure SAWs Monte-Carlo simulations (see Fig. 7 of 
RefQ). 



B. Specific heat 

The specific heat is shown on Fig. 01 for the two boundary conditions : in both cases, the maximum C^} ax scales 
as the size N. For the (bb) case, the temperature T C (N) of the maximum moves towards the critical temperature T c 
as T C (N) — T c ~ l/N, whereas for the (bu) boundary conditions, the temperature of the maximum remains almost 
stable as N varies, as was also found in the pure SAWs Monte-Carlo simulations (see Fig. 8 of Ref [|). 



C. Pairing probabilities along the chain 

On Fig. [51(a), we have plotted the unpairing probability q(a) = 1 — p(a) (see eq. [5jl as a function of the monomer 
index a, at criticality, for both boundary conditions : for the (bb) boundary condition, this unpairing probability is 
flat around the value q c , whereas for the (bu) boundary condition, the unpairing probability varies linearly as 

q(a) ~q c + (l-q c )?L (14) 

This linear behavior actually means that there exists a phase coexistence at T c between the localized phase on (1, ad) 
and the delocalized phase on (ad,N), where the position ad of the interface is uniformly distributed on (1, JV). In 
the SAWs simulations, this coexistence with a uniform position of the interface is seen via the flat histogram of the 
contact density over Monte-Carlo configurations (see Fig. 5 of Ref [l|). For the (bb) boundary condition, there is no 
phase coexistence : the whole sample is in the localized phase. 

For the (bu) boundary condition, the presence of a free end of length I ~ N at criticality means that the choice of 
the function g(a;N) in the free end weight of eq. Q may have some importance for phase coexistence properties at 
criticality. 
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FIG. 4: The specific heat C N of the pure case, for sizes N = 10 5 (O), 2 • 10 5 (□), 4 ■ 10 5 (o), 8 • 10 5 (A). The boundary 
conditions are (a) bound-bound (b) bound-unbound. 
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FIG. 5: (a) The unpairing probability q(a) = (1 — p(a)) of the bound-bound (lower curve) and bound-unbound (upper curve) 
pure case at T c , for TV = 10 5 and ci = (b) The unpairing probability q(a), a > -y, of various (bu) pure cases at T c for 
TV = 10 5 : ci — (A), ci =0.1 (O)i c i — —0.07 (*), and the prescription of reference || (□). Note that the latter is very close 
to the case c\ = studied in this paper. 



A plausible choice, that we have adopted here, is to extend the equation g(a; N) — nv-^w > valid for (N — a) -C N 
to all values of a. We show in Fig 03(b) the unpairing probability q(a) for various values of c\ : c\ ~ 0.1 and 
ci = 0. 

Other choices of the exponent C\ and/or of the function g(a; N) are possible 0>H,[3- We show in particular in Fig 
(b) the unpairing probability q(a) for (i) the value c\ = —0.07, obtained in ref. | 2 1| for a three arm star polymer 
(ii) the prescription g(a; N) = ( N 2 a yi > with c\ ~ 0.1 advocated in ref 

As shown in Fig, 03 we obtain only slight deviations with respect to the straight line behavior of eq. (|14|l obtained 
for our choice g(a; N) = {N \ )C1 , <T = 0. 



D. Conclusion on the pure case 



We have shown how the two boundary conditions (bb) and (bu) lead to different behaviors : at criticality, the (bb) 
case presents a pure localized phase, whereas the (bu) case gives rise to a phase coexistence between a localized phase 
and a delocalized phase, with an interface uniformly distributed over the sample. In this respect, the (bb) case seems 
simpler. However, from the point of view of finite size scaling analysis, the (bu) case has the advantage of allowing a 
direct measure of the critical temperature via the crossing of the contact density. 
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We have also mentioned how the the presents results for the (bu) case can be qualitatively compared with the pure 
SAWs simulations [l|, and we refer to y| for a detailed quantitative comparison between the two models. 
We now turn to the study of the disordered case. 



IV. FINITE SIZE PROPERTIES OF THE DISORDERED CASE 



We consider the case of binary disorder where the binding energies e a are independent random variables: e a = 
Eo = —355 K with probability \ and e Q = E\ = —390 K with probability \. As mentioned in the Introduction, this 
distribution represent AT and GC base pairs, and the values of (eq, £i) are such that the critical temperatures of the 
pure cases are close to the AT and GC experimental ones. 

Disordered averaged quantities A will be denoted by A. 



A. Disordered averaged contact and energy densities 



We have plotted on Fig. the disorder averaged contact density &n(T) for the two boundary conditions : the 
results are very similar to the pure case results of Fig. [21 In particular, for the (bu) case, the crossing of the contact 
density for different values of N gives a direct determination of T c , and shows that the contact density is finite at 
criticality: the transition is therefore first order. We also show in Fig. (b) the contact energy of eq. Q for the (bu) 
case, which is somewhat similar to the corresponding result for SAWs (see Fig 4 of ref The crossing is again a 

clear sign of a discontinuous transition. 
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FIG. 6: The fraction of paired monomers 8 N (T) of the disordered case for sizes N = 10 5 (O), 2 ■ 10 5 (□), 4 • 10 s (A), averaged 
over 10 4 samples. The boundary conditions are (a) bound-bound (b) bound-unbound. The error bars are much smaller than 
the symbols. 



However, in contrast with the pure case, the rescaling (|13Jl shown on Fig. 0(a), is not satisfied at all 

6 N {T)^G{{T-T C )N) (15) 

Since the crossing of the contact density of Fig. (b) implies that the transition is first-order, we are led to the 
conclusion that there is a problem with the finite-size scaling form for the disorder averaged contact density (|15|) . To 
understand why, we have studied the sample-to-sample fluctuations. 



B. Sample to sample fluctuations 



On Fig. |S| (a), we have plotted the contact density for three particular samples, as compared to the average over 
10 4 samples. The sample to sample fluctuations are still strong even for the size N = 10 5 . In particular, at the critical 
temperature T c obtained from the crossing of Fig. El(b), some samples are clearly already 'delocalized', whereas other 
remain localized up to a higher temperature. On Fig. El(b), we have plotted the specific heat for the three same 
samples : the peaks are clearly separated. 
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FIG. 7: (a) Tentative rescaling 9n(T) = 0((T ~T C )N) for the fraction of paired monomers, averaged over 500 samples, for sizes 
N = 10 5 (O), 2 ■ 10 5 (□), 4 • 10 5 (A). The boundary condition is bound-bound (b) The contact energy e N (T) of the disordered 
case for sizes N = 10 5 (O), 2 ■ 10 5 (□), 4 ■ 10 5 (A), averaged over 10 4 samples. The boundary condition is bound-unbound. 
The error bars are much smaller than the symbols. 
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FIG. 8: (a) Comparison of the fraction of paired monomers of three individual disorder samples with the average over 10 4 
samples (O)j f° r s i ze N = 10 5 and bound-bound boundary condition. The thermodynamic critical temperature T c is shown 
(b) The specific heat of the three disorder samples. 



These results are reminiscent of the lack of self-averaging found in second order random critical points |25|, |2jj . In 
these systems, the argument goes as follows [27j (i) off-criticality, the finite correlation length £ allows for a division of 
the sample into independent large sub-samples, which in turn leads to the self-averaging property of thermodynamic 
quantities (here, "independent" means that that the interaction term between these sub-samples is a surface term, 
that can be neglected with respect to volume contributions) (ii) at criticality however, the previous 'subdivision' 
argument breaks down because of the divergence of £. 

In usual first-order transitions, there is no diverging correlation length. The analysis of the disorder effects consists 
in dividing the system into finite sub-systems (or phases), and in taking into account the surface tension between 
phases [l3j]. The first order transition of the PS model is very different for two polymeric reasons. First, there is no 
surface tension between the localized and delocalized phases, as explained in section ITlI CI on the coexistence at the 
pure critical point. Second, there exists here an infinite correlation length, as is very clear in the (bu) case where the 
end segment length diverges as 1/(T C — T) In the (bb) case, this diverging correlation length appears in the loop 
length distribution 

P^(l)^L e -^n withtfT) ~ -i- (16) 
t T^T" J-e — ± 



We are therefore led to the conclusion that the lack of self-averaging we find in the disordered PS model comes from 
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the presence of an infinite correlation length, (even though the transition is first order), that prevents the division of 
the chain into independent sub-chains. 

The strong sample-to sample fluctuations explain why the finite-size scaling analysis Ijl5|l does not work for the 
disordered averaged contact density : each sample (i) has its own pseudo-critical temperature T c (i) where it presents 
a pseudo-first order transition, and the resulting disorder average has no nice finite-size scaling behavior. 



C. Finite-size scaling for the contact density in a given sample 



Since we have obtained that there are strong sample to sample fluctuations that invalidates the usual finite size 
scaling analysis for the disorder av erag ed contact density, we have tried to use more refined finite size scaling theories 
in the presence of disorder |25l l2rl I27L l28q : in these references, it is stressed that the sample to sample fluctuations 
of a given observable are predominantly due to the sample to sample fluctuations of the pseudo-critical temperatures 
T c (i, N). As a consequence, these references express the finite size scaling form of an observable X |27ll28[ as 



X®(T) = N p Qi ((T - T c {i,N))N*) 



(17) 



where T c (i,N) is the pseudo-critical temperature of the sample i, determined for instance as the susceptibility peak 
|27|. or by a minimal distance criterium |28j |. The scaling function Qi a priori also depends on the sample, and we 
refer to |23, UM f° r discussions and examples. 

Here, from the specific heat plotted on Fig |Ht>, we see that beyond the sample dependence of the pseudo-critical 
temperature T c (i), there is also a sample dependent scaling function Qi, in contrast with the case studied in |28| . 
where the curves of various samples were the same up to a translation T — T c (i). 

As a consequence, we have tried to obtain some information from a sample-dependent finite-size scaling analysis 
via the following 'sample replication' procedure : we have compared a given random sample E — (ei, ...,£jv) of size 
N = 10 5 with the sample of size 2N — 2.10 5 obtained by gluing together two copies of the initial sample E, and with 
the sample of size AN = 4.10 5 obtained by gluing together four copies of the initial sample E : the corresponding 
results for the contact density for the two boundary conditions are shown on Fig. [!|] In particular, the crossing 
obtained for the (bu) boundary conditions allows to determine the sample-dependent pseudo-critical temperature 
T c (i). The rescaling of these data with the reduced variable x — (T — T c (i))N is found to be satisfactory in the 
delocalized phase T > T c (i) for the two boundary conditions, see Fig. ^3 We have also plotted on Fig. ^2 the 
associated specific heats : for the two boundary conditions, the maximum value scales as N. 
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FIG. 9: The fraction of paired monomers 8n{T) for a single disorder sample, using the 'sample replication' procedure discussed 
in the text, and for sizes N = 10 s (O)i 2-10 5 (□), 4- 10 5 (A). The boundary conditions are (a) bound-bound (b) bound-unbound. 
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FIG. 10: The scaled fraction of paired monomers 6((T — T C )N)) for the same disorder sample as in Fig El f° r sizes N = 10 5 
(O), 2 • 10 s (□), 4 • 10 5 (A). The boundary conditions are (a) bound-bound (b) bound- unbound. 
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FIG. 11: The specific heat Cjv for the same disorder sample as in Figl^and FigESl for sizes TV = 10 5 (O), 2 ■ 10 5 (□), 4 ■ 10 5 
(A). The boundary conditions are (a) bound-bound (b) bound-unbound. 



D. Pairing probability along the chain in a given sample 

On Fig. I12f a). we have plotted, for a given sample at its pseudo-critical temperature T c (i) determined by Fig. 
(b), the unpairing probability q(a) as a function of the monomer index a for the two boundary conditions. These 
data should be compared with the pure case results of Fig. : beside the disorder fluctuations present in the (bb) 
boundary conditions, the (bu) boundary conditions display an upward trend that shows that there exists some phase 
coexistence between a localized and delocalized phase. To show more clearly this coexistence, we have plotted in Fig. 
I12f b) the difference Aq(a) — qbu{a) — qbb{o) for the same disordered sample, as compared to the linear analog of the 
pure case. 
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FIG. 12: (a) The unpairing probability q(a) for a single disorder sample (i) at its critical temperature T c (i), for size N = 10 5 . 
The boundary conditions are bound-bound (lower curve) and bound-unbound (upper curve) (b) The difference Aq(a) = q^ u (a) — 
qbb{a) between the two previous curves compared with the analog quantity of the pure case at T c (see Figure |S] (a)). 

E. Disorder averaged loop statistics 

We now turn to the probability measure Mn(1) of the loops of length I existing in a sample of size N, as denned 
in eq. ijHJ- In Fig. ED we present the disorder averaged Mjy(0 at criticality, in terms of the rescaled length A = l/N 
for the two boundary conditions (bb) and (bu). In each case, we have added the corresponding pure critical analog 
for comparison. 

In the pure case, the loop length distribution at criticality is given in the thermodynamic limit by the loop weight 
used to defined the PS model (0 P(l) ~ l/l c . Here in the pure finite sample, the exponent c = 2.15 indeed corresponds 
to the slope in the regime I -C N, before the appearance of some curvature near I ~ N. This form for the pure critical 
loop statistics is thus very similar to the corresponding result for pure SAWs 0, 0, In the disordered case, we 
clearly obtain the same slope c = 2.15 in the regime I <C N, and deviations with respect to the pure case only occur in 
the region I ~ N. Our conclusion is thus that in the presence of disorder, the statistics of loops at criticality remains 
P(l) ~ l/l c in the thermodynamic limit as in the pure case. 



20.0 



20.0 




_10.0 I ■ ' ■ ' ■ ' ■ ' -20.0 1 ' 1 ' 1 ' 1 1 1 

-8.0 -6.0 -4.0 -2.0 0.0 -8.0 -6.0 -4.0 -2.0 0.0 



FIG. 13: Log-Log plot of the critical loop distribution Mn(1) in the rescaled variable A = 4j for the disordered case (averaged 
over 10 4 samples) and for the pure case (lower curve). Both systems have size N = 10 5 and the boundary conditions are (a) 
bound-bound (b) bound-unbound. 
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V. CONCLUSION 



In this paper, we have studied the binary disordered Poland-Scheraga model of DNA denaturation, in the regime 
where the pure model displays a first order transition (loop exponent c > 2). Our main conclusion is that the transition 
remains first order in the disordered case, but that disorder averaged observables do not satisfy finite size scaling, as 
a consequence of strong sample to sample fluctuations. 

Let us now consider the relation between the PS model and the Monte Carlo simulations of attracting SAW's. In 
the pure case, the finite-size behavior we have presented for the PS model with (bu) boundary conditions is very 
similar to the corresponding behavior obtained previously in Monte-Carlo simulations of SAW's lj. In the disordered 
case, our result for the crossing of the contact energy (Fig[7|(b)) is somewhat similar to the corresponding results (Fig 
4) of ref. [12j| . The failure of finite size scaling for disordered averaged quantities was also noted in this reference. It 
therefore seems that the PS model with (bu) boundary conditions is an appropriate effective model for SAW's even 
in the presence of disorder. 
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